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A renewal-reward model is developed to predict the 
optimum amount of time that Coast Guard personnel should 
spend investigating a vessel for illicit substances. The 
optimal investigation time is determined with respect to 
three criteria; maximizing the number of arrests, 
maximizing the quantity of drugs confiscated, and 
minimizing the quantity of drugs that escape detection. A 
simulation study indicates that the optimal investigation 
time is very sensitive to underlying distributional 
assumptions. The basic service system model may have wider 
application, i.e., to combat modelling, where it may be 
desirable to investigate a potential target to estimate its 
value before committing limited resources. An adaption of 
the model may also be of help in allocating resources for 
mineral exploration. 
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I. 



INTRODUCTION 



One duty of the United States Coast Guard is the 
interdiction of drug smugglers. When Coast Guard personnel 
board a suspect vessel to search for contraband, a decision 
must be made regarding the amount of time allocated for the 
investigation. If a careful, but lengthy, search is 
conducted, there is a high probability of finding 
contraband if it is present. However, a lengthy 
investigation incurs a penalty. While Coast Guard 
personnel are searching one vessel, other shipping, which 
could have contraband onboard, is passing through the 
patrol area without being examined. At the other extreme, 
if the suspect vessel is examined in a cursory manner, the 
Coast Guard Vessel (CGV) will be able to stop and search 
more shipping but there is a higher probability of not 
discovering the contraband, if present, because of the 
short investigation time. The purpose of this thesis is to 
develop a renewal-reward model to determine the optimum 
investigation time. 
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II. RENEWAL-REWARD MODEL 



Consider the following scenario. All shipping 

encountered by a patrolling CGV is stopped and searched to 
determine if contraband is onboard. After leaving port, a 
CGV searcher a time Si until a vessel is sighted. A 
detected vessel has a probability Pb of being Bad (having 
contraband onboard) and Po = 1 - Pb of being Good (having 
no contraband onboard) . A Bad vessel contains the random 
number J = j(l,2,3,...) units of illicit substance with 

probability bj . The time required to find the first unit 
of the contraband, given that J units are present, is 
modelled here as the minimum of J independent identically 
distributed times having distribution function Fz (t) ; this 
model is illustrative only, and may be altered in various 
realistic directions. If time T is required to locate an 
incriminating unit, then 

P{T>T|J=j| = [l-F2(t)]-l ; j=l,2,3,... 

so, upon removing the condition, 

oo 

1-Ft (t) = P{T>tl = S [l-Fz(t)]Jbj ; j*l,2,3,... (2.1) 

i'l 

Adopt the following decision rule: Establish a 

predetermined investigation time L. When a vessel is 
detected, it is stopped and searched. If no contraband is 
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discovered in the investigation time L, the vessel is 
released and the CGV resumes patrol. If any contraband is 
discovered before the end of the investigation time L, the 
vessel is detained and escorted to port for further search 
and investigation. A time period D is required to escort 
the vessel to base. The patrol cycle ends when the 
escorted vessel arrives at the base. 

The long run average reward per cycle, R, can be 
calculated in the following manner [Ross: pp. 279-294]. 

E[Rc ] 

R = (2.2) 

E[C] 

where: E[Rc] is the expected reward for a patrol 

cycle . 

E[C] is the expected duration of a patrol 
cycle . 

The reward associated with an apprehension can be defined 
in various ways. Three different reward criteria will be 
examined. For each of the three cases, we will determine 
the optimal investigation time L which maximizes the long 
run average cycle reward. 

A. EXPECTED DURATION OF PATROL CYCLE: E[C] 

We must distinguish between the following two cases for 
the investigation time I: 

Ig - an investigation time that results in the release of 
the vessel. 

Ib - an investigation time that results in the detention 
of the vessel. 
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Where 



I 



Ig = L with probability Pg + Pb Ft (L) = a 

Ib = T with probability Pb Ft (L) = oT 



The length of a CGV patrol cycle can be represented as: 



where C’ has the same unconditional distribution as C, for 
if the first investigation results in a release of the 
vessel being searched, the process re-starts (regenerates) . 
The expected length of a patrol cycle can be expressed 

as : 



C 



Si + Ib + D With probability a(L) 

Si + Ig + C With probability a(L) 



(2.3) 



E[C] = E[S] + (E[Ibl + E[D])a + (L + E[C'])a 



or 



E[C1 



E[S] + E[Dla + E[Ib1a + La 
1 - a 



(2.4) 



which can be written as 



E[C] = E[D] + E[Ibl + 




(2.5) 



Solving for the term E[Ibl yields 




( 2 . 6 ) 



E[Ib 1 



Ft (L) 



Integration by parts shows that 




tF(dt) 



(2.7) 
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Hence 



Ft (t)dt - LFt (L) 

E[Ib] = (2.8) 

Ft (L) 



Therefore 



E[C] 



E[S] + E[D]PbFT (L) 



+ 




Ft (t)dt + LPfl 



PbFT (L) 



(2.9) 



which reduces to 



E[C] 



E[D] 



+ ECS] + 
PbFT (L) 



Ft (t)dt ^ 
Ft (L) 



Pg L 

PbFT (L) 



( 2 . 10 ) 



B. REWARD CRITERIA 1 

Suppose the CGV is rewarded for making an arrest. 
Since only one arrest is made per cycle. 



E[Rc] = 1 



( 2 . 11 ) 



and from equation 2.2 



R 



_l 

E[C] 



( 2 . 12 ) 



5 



Consider the following example 



F2 (t) 


= 1 - e-^" 


00 

o 

• 

o 

II 


S(t) 


= 1 - e“^^ ^ 


, Xi =1.1 


D(t) 


= 1 - e-^^t 


, X 2 =0.3 


bj = 


(1 - p)p^"^ 


II 

0 

• 

00 

(- 1 . 


Pb = 


0.2 





From equation 2.1 

~ e-ut (1-p) 

Ft (t) = S (e-i'MJ (l-3)pJ-i = (2.13) 

i=l 1 - Pe-«t 

and 

1 - e-M‘ 

Ft (t) = 1 - Ft (t) = (2.14) 

1 - Pe-ut 



From equation 2.10, the expected duration of the patrol 
cycle is 




+ 



Pg L 

P. (l— ) 



(2.15) 



and from equation 2.12, the long run average reward per 
cycle can be calculated. 
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R = 



1-e-uL 






Xi X 2 Pb 



(2.16) 



l-e-ML \ /l-pe-ML 

Pb I )Xi + Xz + Xi Xz Pb I |Ln I )+ PgLXiXz 



l-Pe-jiL/ 



PU y y 1-p 



C. REWARD CRITERIA 2 

The reward assigned to the CGV is the amount of drugs 
confiscated. 



E[Rc] = E[J|T ^ L] 



(2.17) 



The probability that a vessel being searched contains J 
j units of drugs and that the drugs are discovered is 

P{ (T^L)^(J=j) ) = [1 ■ (1-Fz (L) )-ilbj (2.18) 



SO 



= - [1 - 

3=1 



P{T^L1«Ft (L) = Z: 1 - (l-Fz(L))J bj (2.19) 



and thus 



1^1 - (1-Fz (L) )J jbj 



P{ (J=j) I (TSL) I = 



( 2 . 20 ) 



Ft (L) 



and 



E[J|T^L] = 



Z: j [1 - (1-Fz (L) ) j]bj 
j=l *- 

Ft (L) 



( 2 . 21 ) 
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using the parameters from the Example in Model 1, 



S j 1 - (e-'“- (l-e)p-' - ^ 

j.l L J 

E[J|F^L] = (2.22) 

(1-e-KL y ) 

which can be evaluated as in the previous example as 



(1-p) 



E[JIT^L] = 



(^ / (1-P)2^ /(l-pe-***- ] 

(1-P)) " /d-pe-i-'- )) ] 



(2.23) 



using equation 2.2, the long run average reward is 



^ _ E[J T^L3 
^ " E[C] 



or 



R = 



^(l-p)V \ (i-pe-uL)*/ J 

(— ) - 1 

\ (1-P) / V (I-Pe-VL ) J J 



(2.24) 



1 

X 2 



Pb 



1_ 

Xi 

l_e-ML 

l-pe-uL 









1-P 



PcL 



l_e-«L 



Pb 



1-pe- 



U L 



, 1-e-ML V 

\ I-Pe-UL/ 



D. REWARD CRITERIA 3 

The reward examined in model three is the difference 
between the quantity of drugs confiscated and the quantity 
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of drugs that escape detection. The expected reward in a 
cycle is 



Rc = E[J|T5L] - E[M] 



where M is defined as the quantity of drugs that escape 
detection in a cycle. 



Jib = Drug quantity passing through patrol area while CGV 
searches a Bad ship 

Jd = Drug quantity passing through patrol area while CGV 
escorts a Bad ship to shore 

Jib = Drug quantity passing through patrol area while CGV 
searches a Good ship 

Jo = Drug quantity onboard a ship classified as Good 
that is actually Bad 



If conditional expectations are taken assuming that the Bad 
vessels pass through the region according to a Poisson 
process with rate Xi Pb , then 



E[M] =* XiPbE[J] (E[Ib]+E[D])PbFT (L) + (XiPbE[J]E[Ifl]+E[M]) 



M 



Jib + Jd 
Jig + Jo + M 



With probability PbFr(L) 

With probability Po + PbFt (L) 



where : 



(Po+PbFr (L) )+E[Jo iPbFx (L) 



( 2 . 25 ) 
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or 



E[M] =. 



Xi PbE[J] Pb 



i: 



Ft (Ddt-LFb 



^ +E[D] Pb 



Ft (L)+L Pg +PbFT 



Pb Ft (L) 



E[Jo ]PbFT (L) 



PbFT (L) 



where E[Jo] is defined as 



(2.26) 



E[J|T>L] = 



S j (l-Fz (L) )J bj 

j=l 

fT (L) 



Using the distributions defined in the Example for Model 1 



E[M] = 



PbXi 



r 1-p /l-pe-i^'-v i-e-w*- 

Ln 1+ 



(1-p) 



1-e-uL 



LPg 






1-3 



I-Pe-ML / 

(1-P)e-M‘- 



(1-Pe-“L) (i-e-ML) 



A2(l-3e->“-) Pb 



(2.27) 



and the long run average reward can be expressed as 



R = Eguation [2.21] - Equation [2.27] (2.28) 

Equation [2.15] 
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III. RESULTS 



The optimum investigation time was determined for the 
above example for each reward criteria using two methods. 
The first method consisted of writing a Monte Carlo 
simulation of the patrol cycle using the distributions of 
the example and running the simulation using various values 
of the investigation time L. Statistics were gathered 
during the simulation allowing the calculation of the long 
run average reward R. The investigation time was varied 
from zero to six hours in increments of two tenths of an 
hour. The optimum investigation time was determined by 
graphing the long run average reward as a function of 
investigation time and finding the value of L which 
maximized R. Each cycle was replicated 20,000 times. A 
detailed discussion of the simulation can be found in 
Appendix A. 

THe second method consisted of writing a computer 
program for the three equations representing the long run 
average reward for the three reward criteria for the 
example and solving the equations for various values of the 
investigation time L. As in the case of the simulation, 
the value of L was varied from zero to six hours in 
increments of two tenths of an hour. The optimum 
investigation time was again determined by graphing the 
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long run average reward as a function of investigation time 
and finding the value of L which maximized R. 

Figures 1,2, and 3 contain the results of both the 
simulation and numerical solution using the distributions 
of the examples for reward criteria 1,2, and 3 
respectively. The solid line represents the analytical 
solution and the circles are the simulation results. The 
maximum long run average reward using criteria 1 can be 
achieved by using an investigation time between 1.4 and 2.2 
hours. The maximum long run average reward using criteria 
2 and 3 can be achieved using an investigation time between 
1.2 and 1.4 hours. It is interesting to note that, using 
the distributions presented in these examples, and 
investigation time exists that maximizes all three reward 
criteria simultaneously. This occurs at 1.4 hours for the 
input distributions. We do not know that this state of 
affairs will persist. 
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Figure 1: Reward-Criteria 1 Results 
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Figure 2: Reward Criteria 2 Results 
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( E(J) - E(M) ) / E(C) 



O 




INVESTIGATION TIME 

Figure 3: Reward Criteria 3 Results 
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IV. SENSITIVITY ANALYSIS 



In order to check the robustness of the model with 
respect to the underlying assumptions, the simulation was 
run exactly as before with the exception that Fz (t) is now 
assumed to have a lognormal distribution instead of an 
exponential distribution. Each cycle was replicated 20,000 
times. Three separate cases are examined. In case A, the 
distribution function Fz (t) has the same mean and variance 
as the exponential distribution used in the previous 
examples. In case B, the distribution function Fz (t) has 
the same mean but twice the variance as the exponential 
distribution used in the previous examples. Finally, in 
case C, the distribution function Fz (t) has the same mean 
but four times the variance as the exponential distribution 
used in the previous examples. 

Figures 4,5, and 6 contain a comparison between the 
exponential base case and the lognormal cases mentioned 
above for reward criteria 1,2, and 3 respectively. It is 
readily apparent that the results for the exponential base 
case and lognormal case A, with the same mean and variance, 
are different. Furthermore, within the lognormal family of 
curves, it can be seen by examining case B and C that the 
results vary significantly as the variance increases. 
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Figure 4; Reward Criteria 1 Sensitivity Analysis 
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Figure 5: Reward Criteria 2 Sensitivity Analysis 
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Figure 6: Reward Criteria 3 Sensitivity Analysis 
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In order 


to better 


understand 


this behavior, it is 


helpful 


to compare the 


quantiles of 


T, the time 


required to 


locate 


an 


incriminating unit given that J 


units are 


present , 


for 


the four distributions. 






Quantile 




Base Case 


Case A 


Case B 


Case C 


0.10 




0.242 


1.604 


0.908 


0.431 


0.25 




0.852 


2.502 


1.467 


1.681 


0.50 




2.465 


3.945 


2.856 


1.681 


0.75 




6.504 


6.858 


5.360 


3.909 


0.90 




13.488 


12.155 


10.333 


11.245 



It can be seen that some distributions, in particular 
the lognormal, cases a, b, and c, have a far greater 
likelihood of producing large values of T relative to the 
others, such as the exponential distribution used in the 
examples. Since the time required to find the first unit 
of drugs greatly influences the cycle length, distributions 
generating larger values of T will produce significantly 
different results. To further verify the accuracy of the 
distributions produced by the simulation and the quantiles 
listed above. Appendix B contains the calculated quantiles 
for lognormal case A produced using analytical distribution 
theory . 
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V. 



CONCLUSIONS 



It has been demonstrated that a renewal-reward approach 
to modelling the Coast Guard drug interdiction process is 
feasible and that it ds possible to determine an optimal 
investigation time. By considering several different 
reward criteria in the model, it is possible to quantify 
the costs with respect to the rest of the reward criteria 
when one criteria is selected and used to arrive at an 
optimum investigation time. 

Of the two methods presented to obtain the optimum 
investigation time, the simulation approach is the most 
flexible. The example distributions and assumptions used 
in this thesis were selected because they allowed an 
analytical solution to be found. This allowed a comparison 
of results between the simulation and numerical solutions, 
thus verifying the simulation. As the model assumptions 
are changed and different probability distributions 
incorporated, the simulation can be easily modified to 
reflect these changes whereas an analytical solution may no 
longer be possible. 

As the model is currently formulated, its usefulness is 
questionable due to the sensitivity to the underlying 
distribution of the time required to find drugs onboard a 
vessel given drugs are present. As demonstrated in the 
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sensitivity analysis even within a given family of 
distributions, the model is sensitive to the distribution 
parameters. In order for the model to produce realistic 
results, rigorous data analysis must be conducted to 
properly identify this distribution. 

The basic service system model may have wider 
application, e.g., to combat modelling, where it may be 
desirable to investigate a potential target to estimate its 
value before committing limited resources. An adaption of 
the model may also be of help in allocating resources for 
mineral exploration. 
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APPENDIX A 



The simulation is written in Fortran 77. It consists 
of a main program and seven subroutines. Uniform [0,1] 
random numbers are provided by calling procedure GGUBFS 
located in the IBM IMSL single precision library. All real 
variables are computed using double precision to minimize 
rounding error. Output is directed to three separate 
units. Unit 2 contains the calculated rewards for each 
increment of investigation time using reward criteria 1, 2, 
and 3. Unit 3 contains the cycle length, drug quantity 
confiscated, and quantity of drugs missed for each 
increment of investigation time. Unit 4 contains detailed 
information regarding the cycle for each increment of 
investigation time. This information includes the total 
number of ships searched, the number of ships that are 
good, bad, bad but declared good, and bad identified as 
bad. Also included in unit 4 output is the time required 
to discover the drugs on a ship declared bad, the quantity 
of drugs onboard ships passing through the area while the 
CGV is searching a ship, the quantity of drugs missed due 
to short investigation time, and the quantity of drugs 
confiscated on a bad ship. 

The main program controls the starting and final 
investigation time, investigation time increment width, and 
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the number of cycles per investigation time to be 
simulated. Cycle averages are computed and output directed 
to the three units discussed above. 

Subroutine CYCLE simulates one patrol cycle and records 
all the relevant statistics during the cycle. Subroutine 
STIME generates the random search times required to find a 
ship from an exponential distribution. Subroutine DTIME 
generates the random times required to escort a ship back 
to base from an exponential distribution. Subroutine 
DRGQTY generates the random quantity of drugs on a BAD ship 
from a geometric distribution. Subroutine CLASS determines 
the classification of a vessel based on the deterministic 
value Pb . Subroutine OPCOST determines the quantity of 
drugs missed onboard other ships passing through the area 
while the CGV is investigating the current ship. 
Subroutine RTIME generates the random times required to 
find the first unit of drugs given that j units are 
present . 

Four versions of this subroutine are listed; one for 
the exponential distribution and the other three for the 
lognormal case A, B, and C distributions. 
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FILE: DPSIM 



FORTRAN A1 



THIS PROGRAM SIMULATES PATROL CYCLES FOR VARIOUS LENGTHS OF 
INVESTIGATION TIMES AND CALCULATES THE VALUES FOR THE THREE 
REHARD CRITERIA. FUNCTION GGUBFS FROM THE IMSLSP LIBRARY IS 
USED TO GENERATE UNIFORM (0,1) NUMBERS. 

INTEGER I,K,MAXT,NCYCLE,Z 

REALMS L,C, J,M,INC,EC(500),EJ(500),EM(500),ITIME(500),TOEND, 

& R1(500),R2(500),R3(500),TS,GS,BS,BSG,BSB,EIB,OP,MISS,EIG, 

& CATCH, ES, ED, BJ,DSEED 
DSEED = 995317. DO 
MAXT = A 
NCYCLE = 100 
INC = 0.20 
L = 0 

TOEND = MAXT/INC 
DO 200 K = 1, TOEND 

L = L + INC 

C = 0 

J = 0 

M = 0 

PRINTX, 'COMPUTING L = *,L 
DO 100 I = 1, NCYCLE 

CALL CYCLE(DSEED,L,C, J,M,TS,GS,BS,BSG,BSB,EIB,OP,MISS, 
& CATCH, ES, ED, BJ,EIG) 

100 CONTINUE 

Z=NCYCLE 
NRITE(A,x) L 

NRITE(A,523) TS/Z, GS/Z, BS/Z, BSG/Z, BSB/Z 
NRITE(A,52A) EIB/Z, B J/Z, OP/Z, MISS/Z, CATCH/Z 
WRITE(A,525) EIG/Z, ES/Z, ED/Z 
WRITE(A,x) 

523 F0RMAT(1X,5F9.A) 

52A F0RMAT(1X,5F9.A) 

525 F0RMAT(1X,3F9.4) 

TS =0 
BS =0 
BSG=0 
BSB=0 
GS =0 
EIB = 0 
MISS=0 
CATCH=0 
ES =0 
ED =0 
BJ =0 
OP =0 
EIG=0 

ITIME(K) = L 
EC(K) = C/NCYCLE 
EJ(K) = J/NCYCLE 
EM(K) = M/ NCYCLE 
Rl(K) = l/ECCK) 

R2(K) = EJ(K)/EC(K) 

R3(K) = (EJ(K)-EM(K))/EC(K) 

NRITE(2,19) ITIME(K),R1(K),R2(IC),R3(K) 

NRITE(3,19) ITIME(K),EC(K),EJ(K),EM(K) 

19 FORMAT (1X,F6.2,3F9.3) 

200 CONTINUE 

STOP 'PROGRAM COMPLETE' 

END 



SUBROUTINE CYCLE(DX,LSTAR,CSTAR, JSTAR,MSTAR,TS,GS,BS,BSG,BSB, 

& EIB,OP,MISS,CATCH,ES,ED,BJ,EIG). 

REALMS LSTAR,CSTAR,JSTAR,MSTAR,LAMDA1,LAMDA2,MU,BETA,S,R,QTY, 

S D, LEAK, TIME, PBAD,TS,GS,BSG,BSB,EIB, OP, MISS, ES, ED, BJ,BS, CATCH, 
a EIG,DX 
INTEGER NC 
CHARACTERX4 TYPE 
LAMDAl =1.10 
LAMDA2 =0.3 
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MU = 0.08 
BETA =0.8 
PBAD = 0.2 

20 TS=TS+1 

CALL STIME(DX,LAMDA1,S) 

ES=ES+S 

CALL CLASSCDX, PBAD, TYPE) 

IFC TYPE .EQ. 'GOOD*) THEN 
EIG = EIG+LSTAR 
GS=GS+1 
QTY=0 

TIME = LSTAR 

CSTAR = CSTAR + TIME + S 

CALL OPCOST(DX,LAMDAl,TIME, BETA, PBAD, LEAK) 

MSTAR = MSTAR + LEAK 

OP=OP+LEAK 

GOTO 20 

ELSE IFC TYPE .EQ. *BAD») THEN 
BS=BS+1 

CALL DRGQTY(DX,BETA,QTY) 

BJ=BJ+QTY 

CALL RTIME(DX,MU,QTY,R) 

IFCR .GT. LSTAR) THEN 
EIG=EIG+LSTAR 
BSG=BSG+1 
TIME = LSTAR 

CSTAR = CSTAR + TIME + S 

CALL OPCOST(DX,LAMDAl,TIME, BETA, PBAD, LEAK) 

MSTAR = MSTAR + LEAK + QTY 

OP=OP+LEAK 

MISS=MISS+QTY 

GOTO 20 

ELSEIFCR .LE. LSTAR) THEN 
EIB=EIB+R 
BSB=BSB+1 

CALL DTIME(DX,LAMDA2,D) 

ED=ED+D 

TIME = R + D 

CSTAR = CSTAR + TIME + S 

JSTAR = JSTAR + QTY 

CALL OPCOST(DX,LAMDAl, TIME, BETA, PBAD, LEAK) 

MSTAR = MSTAR + LEAK 

OP=OP+LEAK 

CATCH=CATCH+QTY 

GOTO 21 

ELSE 

STOP 'ERROR 1* 

ENDIF 

ELSE 

STOP 'ERROR 2* 

ENDIF 

21 RETURN 
END 



SUBROUTINE STIMEC DRSEED, INPUT, OUTPUT) 
REALX8 INPUT, OUTPUT, RV, DRSEED 
RV=GGUBF3(DRSEED) 

OUTPUT = LOGCRV)/(-INPUT) 

RETURN 

END 



SUBROUTINE DTIMEC DRSEED, INPUT, OUTPUT) 
REALX8 INPUT, OUTPUT, RV, DRSEED 
RV = GGUBFS(DRSEED) 

OUTPUT = LOG(RV)/(-INPUT) 

RETURN 

END 
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SUBROUTINE DRGOTYC DRSEED, INPUT, OUTPUT) 
REAL3^8 DRSEED, INPUT, OUTPUT, RV 
RV = GGUBFS(DRSEED) 

OUTPUT = AINT(LOG(RV)/(LOG(INPUT)))+l 

RETURN 

END 



SUBROUTINE CL ASS ( DRSEED, INPUT, OUTPUT) 
REALX8 DRSEED, INPUT, RV 
CHARACTER^^A OUTPUT 
RV = GGUBFS(DRSEED) 

IFCRV .LE. INPUT) THEN 
OUTPUT = 'BAD' 

ELSE 

OUTPUT = 'GOOD' 

ENDIF 

RETURN 

END 



SUBROUTINE OPCOSTC DRSEED, INPUTl , INPUT2, PI , P2, OUTPUT) 
REAL5^8 DRSEED, INPUTl, INPUT2,P1,P2, OUTPUT, STOR,NSHIP 
OUTPUT = 0 
JO = 1/(1-P1) 

NSHIP = INPUT25^INPUT1 
OUTPUT = (P2^NSHIP))CJ0 
RETURN 
END 



SUBROUTINE RTIMEC DRSEED, INPUTl , INPUT2, OUTPUT ) 

C EXPONENTIAL DISTRIBUTION - BASE MEAN S VARIANCE 
REALMS DRSEED, INPUTl, INPUT2, OUTPUT, RV 
RV = GGUBFS(DRSEED) 

OUTPUT = LOG(RV)/(-INPUTl9«INPUT2) 

RETURN 

END 



SUBROUTINE RTIMEl ( DRSEED, INPUTl , INPUT2, OUTPUT) 

C LOGNORMAL DISTRIBUTION - BASE MEAN S VARIANCE 

REALX8 INPUTl, INPUT2, OUTPUT, STDV,AVG, PI, LOW, U1,U2,GG, 
& DRSEED 

INTEGER II, FINISH 
GG = LOGC2.0) 

LOW = 999999 
FINISH = INTCINPUT2) 

STDEV = SQRT(GG) 

AVG = -(LOG(INPUTl))- O.S^GG 
PI = 3.1A1592654 
DO 662 II = 1, FINISH 
U1 = GGUBFS(DRSEED) 

U2 = GGUBFS(DRSEED) 

RV = S0RT(-2^L0G(U1))^C0S(2^PD€U2)^STDEV+AVG 
RVT = EXP(RV) 

IFCRVT .LT. LOW) THEN 
LOW = RVT 

ENDIF 

662 CONTINUE 

OUTPUT = LOW 

RETURN 

END 

SUBROUTINE RTIME2C DRSEED, INPUTl , INPUT2, OUTPUT) 

C LOGNROMAL DISTRIBUTION - BASE MEAN & 2^VARIANCE 

REALMS INPUTl, INPUT2, OUTPUT, STDV, AVG, PI, LOW, U1,U2,GG, 
& DRSEED- 
INTEGER II, FINISH 
GG = LOGC3.0) 

LOW = 999999 
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FINISH = INTCINPUT2) 

STDEV = SQRT(GG) 

AVG = -(LOG(INPUTl))- O.S^GG 
PI = 3.1A159265A 
DO 662 II = 1, FINISH 
U1 = GGUBFS(DRSEED) 

U2 = GGUBFS(DRSEED) 

RV = SQRT(-25^L0G(U1))5«C0S(23<PI5(U2)X$TDEV+AVG 
RVT = EXP(RV) 

IFCRVT .LT. LOW) THEN 
LOW = RVT 

ENDIF 

662 CONTINUE 

OUTPUT = LOW 

RETURN 

END 

SUBROUTINE RTIME3C DRSEED, INPUTl , INPUT2, OUTPUT) 

C LOGNROMAL DISTRIBUTION - BASE MEAN & AXVARIANCE 

REALMS INPUT1,INPUT2,0UTPUT,STDV,AVG>PI,L0W,U1,U2,GG, 
& DRSEED 

INTEGER II, FINISH 
GG = LOGC5.0) 

LOW = 999999 
FINISH = INTCINPUT2) 

STDEV = SQRT(GG) 

AVG = -(LOG(INPUTl))- 0.5XGG 
PI = 3.1A159265A 
DO 662 II = 1, FINISH 
U1 = GGUBFS(DRSEED) 

U2 = GGUBFS(DRSEED) 

RV = SQRT(-2)^L0G(U1))XC0S(2XPIXU2)5(STDEV+AVG 
RVT = EXP(RV) 

IFCRVT .LT. LOW) THEN 
LOW = RVT 

ENDIF 

662 CONTINUE 

OUTPUT = LOW 

RETURN 

END 
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THIS PROGRAM COMPUTES THE VALUES OF THE THREE REWARD CRITERIA 
FOR VARIOUS VALUES OF L FOR COMPARISON WITH SIMULATION RESULTS 

REAL PB,LAMDA1,LAMDA2,BETA,MU,INC,ES,ED,FT,FTBAR,A,ABAR, 
a L,TOEND,EC(100),EJ(100),EM(100),R1(100),R2(100),R3(100). 
a ITIME(100),PG,B,BB,E,C1.C2,C3,CA,M1,M2,TEMP,J1,RATE 
INTEGER MAXT,K 
MAXT = 6 
INC = 0.20 
PB = .20 
PG = 1-PB 
LAMDAl =1.10 
LAMDA2 =0.3 
E5=1/LAMDA1 
ED=1/LAMDA2 
BETA =0.8 
MU = 0.08 
B=BETA 
BB=1-B 

RATE = (LAMDA1XPB)/BB 
TOEND = MAXT/INC 
L=0 

DO 100 K = 1. TOEND 
L=L+INC 
E = EXP(-MU3^L) 

TEMP = 1-(B3^E) 

FTBAR = (BBXE)/TEMP 
FT=1-FTBAR 
ABAR=PB?^FT 
A=1-ABAR 
Cl= ES/ABAR 
C2=ED 

C3 = (PB/ABAR)?^((BB/(B5^MU))5(LOG(TEMP/BB)-(L5^ FTBAR)) 

CA = (L3^A)/ABAR 

J1 = (BB/FT)3<( (l/(BB)fBB))-(E/(TEMP3^TEMP))) 

M1 = (C3+ED))^RATE+((L5(A9^RATE)/ABAR) 

M2 = BB5^(E/(TEMPXTEMP))X(1/ABAR)XPB 
ITIME(K) = L 
EC(K) = C1+C2+C3+C4 
EJ(K) = J1 
EM(K) = M1+M2 
Rl(K) = l/ECCK) 

R2(K) = EJ(K)/EC(K) 

R3(K) = (EJ(K)-EM(K))/EC(K) 

WRITE(2,132) ITIME(K),R1(K).R2(K),R3(K) 

WRITE(3,132) ITIME(K),EC(K)^EJ(K),EM(K) 

WRITE(A,3«) L 

WRITE(A,133) C1,C2,C3.C4 
WRITE(A,13A) J1 
WRITE(^,135) Ml, M2 

132 FORMAT ( IX, F6 . 2, 3F9 . 3) 

133 FORMAT (1X,AF12.5) 

134 FORMAT (1X,F12.5) 

135 FORMAT (1X,2F12.5) 

WRITE(4,X) 

100 CONTINUE 

STOP 'PROGRAM EXECUTION COMPLETE' 

END 
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APPENDIX B 



Recall that the probability distribution for 
quantity of drugs on a vessel given that it is bad is: 



the 



Let 



so 



bj 



= (1 - 
P{T>t t 

Ft (t) 



3 ) PJ - 1 ; 3 = 



s Ft (t) = S 

j=l 



(l-B)Fz (t) 
1-BFz (t) 



0.8, j = 1,2,3,... 
(l-Fz (t) )Jbj 



and 



P{T^tl 



(1-B) Fz (t) ^ 1-Fz_(t) 

l-BFT{t) 1-BFz (t) 



now let 



1-Fz (tpj ^ 
1-BFz <tp ) 



which can be written as 



and 



Fz (tp ) 



1-P 

I-Pp 



Fz (tp ) = 1 - Fz (tp ) = 
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Assuming Fz (t) is lognormal with mean = 2.1792 and standard 



deviation = 0.8326 as in distribution Case A 




In tp - u 



a 



P(l-B) 

1-Bp 



therefore 



tp — EXP p + CZp ( 1 - D > 

1 - 0 p 



] 



where Z p ( i - p > is taken from tables of the Standard Normal 

1 - P p 

distribution. 

Using the equation for tp , the following table of quantiles 
was generated: 



It can be seen that these quantile values agree closely 
with those presented in the sensitivity analysis. 



P 



tp 



0.10 

0.25 

0.50 

0.75 

0.90 



1.65 

2.46 

3.95 

6.78 

11.99 
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